home *** CD-ROM | disk | FTP | other *** search
/ CD Actual 1 / PC Actual CD 01.iso / share / dos / graficos / plydat14.arj / LIB.C < prev    next >
Encoding:
C/C++ Source or Header  |  1992-03-29  |  14.6 KB  |  555 lines
  1. /*
  2.  * lib.c - a library of vector operations, a random number generator, and
  3.  *     object output routines.
  4.  *
  5.  * Version:  2.2 (11/17/87)
  6.  * Author:  Eric Haines, 3D/Eye, Inc.
  7.  */
  8.  
  9. #include <stdio.h>
  10. #include <stdlib.h>
  11. #include <math.h>
  12. #include "def.h"
  13. #include "lib.h"
  14.  
  15. static texture_count = 0;
  16. /*
  17.  * Normalize the vector (X,Y,Z) so that X*X + Y*Y + Z*Z = 1.
  18.  *
  19.  * The normalization divisor is returned.  If the divisor is zero, no
  20.  * normalization occurs.
  21.  *
  22.  */
  23. double
  24. lib_normalize_coord3(COORD4 *cvec)
  25. {
  26.     double divisor;
  27.  
  28.  
  29.     divisor = sqrt( (double)DOT_PRODUCT( (*cvec), (*cvec) ) ) ;
  30.  
  31.     if ( divisor != 0.0 ) {
  32.         cvec->x /= divisor;
  33.         cvec->y /= divisor;
  34.         cvec->z /= divisor;
  35.     }
  36.  
  37.     return( divisor );
  38. }
  39.  
  40.  
  41. /*
  42.  * Set all matrix elements to zero.
  43.  */
  44. void
  45. lib_zero_matrix(MATRIX mx )
  46. {
  47.     int   i, j ;
  48.     for ( i = 0 ; i < 4 ; ++i ) {
  49.         for ( j = 0 ; j < 4 ; ++j ) {
  50.             mx[i][j] = 0.0 ;
  51.         }
  52.     }
  53. }
  54.  
  55.  
  56. /*
  57.  * Create identity matrix.
  58.  */
  59. void
  60. lib_create_identity_matrix(MATRIX mx)
  61. {
  62.     int i;
  63.  
  64.  
  65.     lib_zero_matrix( mx ) ;
  66.     for ( i = 0 ; i < 4 ; ++i ) {
  67.         mx[i][i] = 1.0 ;
  68.     }
  69. }
  70.  
  71. /* Create translation matrix */
  72. void
  73. lib_create_translate_matrix(MATRIX mx, COORD4 *vec)
  74. {
  75.     int i;
  76.     lib_create_identity_matrix(mx);
  77.     mx[3][0] = vec->x;
  78.     mx[3][1] = vec->y;
  79.     mx[3][2] = vec->z;
  80. }
  81.  
  82. /* Create translation matrix */
  83. void
  84. lib_create_scale_matrix(MATRIX mx, COORD4 *vec)
  85. {
  86.     int i;
  87.     lib_zero_matrix(mx);
  88.     mx[0][0] = vec->x;
  89.     mx[1][1] = vec->y;
  90.     mx[2][2] = vec->z;
  91. }
  92.  
  93. /*
  94.  * Create a rotation matrix along the given axis by the given angle in radians.
  95.  */
  96. void
  97. lib_create_rotate_matrix(MATRIX mx, int axis, double angle)
  98. {
  99.     double  cosine ;
  100.     double  sine ;
  101.  
  102.  
  103.     lib_zero_matrix( mx ) ;
  104.  
  105.     cosine = cos( (double)angle ) ;
  106.     sine = sin( (double)angle ) ;
  107.  
  108.     switch ( axis ) {
  109.         case X_AXIS:
  110.             mx[0][0] = 1.0 ;
  111.             mx[1][1] = cosine ;
  112.             mx[2][2] = cosine ;
  113.             mx[1][2] = sine ;
  114.             mx[2][1] = -sine ;
  115.             break ;
  116.         case Y_AXIS:
  117.             mx[1][1] = 1.0 ;
  118.             mx[0][0] = cosine ;
  119.             mx[2][2] = cosine ;
  120.             mx[2][0] = sine ;
  121.             mx[0][2] = -sine ;
  122.             break ;
  123.         case Z_AXIS:
  124.             mx[2][2] = 1.0 ;
  125.             mx[0][0] = cosine ;
  126.             mx[1][1] = cosine ;
  127.             mx[0][1] = sine ;
  128.             mx[1][0] = -sine ;
  129.             break ;
  130.     }
  131.     mx[3][3] = 1.0 ;
  132. }
  133.  
  134.  
  135. /*
  136.  * Create a rotation matrix along the given axis by the given angle in radians.
  137.  * The axis is a set of direction cosines.
  138.  */
  139. void
  140. lib_create_axis_rotate_matrix(MATRIX mx, COORD4 *rvec, double angle)
  141. {
  142.     COORD4  axis ;
  143.     double  cosine ;
  144.     double  one_minus_cosine ;
  145.     double  sine ;
  146.  
  147.  
  148.     lib_zero_matrix( mx ) ;
  149.  
  150.     COPY_COORD( axis, (*rvec) ) ;
  151.  
  152.     cosine = cos( (double)angle ) ;
  153.     sine = sin( (double)angle ) ;
  154.     one_minus_cosine = 1.0 - cosine ;
  155.  
  156.     mx[0][0] = SQR(axis.x) + (1.0 - SQR(axis.x)) * cosine ;
  157.     mx[0][1] = axis.x * axis.y * one_minus_cosine + axis.z * sine ;
  158.     mx[0][2] = axis.x * axis.z * one_minus_cosine - axis.y * sine ;
  159.  
  160.     mx[1][0] = axis.x * axis.y * one_minus_cosine - axis.z * sine ;
  161.     mx[1][1] = SQR(axis.y) + (1.0 - SQR(axis.y)) * cosine ;
  162.     mx[1][2] = axis.y * axis.z * one_minus_cosine + axis.x * sine ;
  163.  
  164.     mx[2][0] = axis.x * axis.z * one_minus_cosine + axis.y * sine ;
  165.     mx[2][1] = axis.y * axis.z * one_minus_cosine - axis.x * sine ;
  166.     mx[2][2] = SQR(axis.z) + (1.0 - SQR(axis.z)) * cosine ;
  167.  
  168.     mx[3][3] = 1.0 ;
  169. }
  170.  
  171. /* Given a point and a direction, find the transform that brings a point
  172.    in a canonical coordinate system into a coordinate system defined by
  173.    that position and direction. */
  174. void
  175. lib_create_canonical_matrix(MATRIX trans, COORD4 *origin, COORD4 *up)
  176. {
  177.    MATRIX trans1, trans2;
  178.    COORD4 tmpv;
  179.  
  180.    /* Translate "origin" to <0, 0, 0> */
  181.    lib_create_translate_matrix(trans1, origin);
  182.  
  183.    /* Determine the axis to rotate about */
  184.    if (fabs(up->z) == 1.0)
  185.       SET_COORD4(tmpv, 1.0, 0.0, 0.0, 1.0)
  186.    else
  187.       SET_COORD4(tmpv, -up->y, up->x, 0.0, 1.0)
  188.    lib_normalize_coord3(&tmpv);
  189.    lib_create_axis_rotate_matrix(trans2, &tmpv, acos(up->z));
  190.    lib_matrix_multiply(trans, trans2, trans1);
  191. }
  192.  
  193.  
  194. /*
  195.  * Multiply a 4 element vector by a matrix.
  196.  */
  197. void
  198. lib_transform_coord(COORD4 *vres, COORD4 *vec, MATRIX mx)
  199. {
  200.     vres->x =
  201.         vec->x*mx[0][0] + vec->y*mx[1][0] + vec->z*mx[2][0] + vec->w*mx[3][0] ;
  202.     vres->y =
  203.         vec->x*mx[0][1] + vec->y*mx[1][1] + vec->z*mx[2][1] + vec->w*mx[3][1] ;
  204.     vres->z =
  205.         vec->x*mx[0][2] + vec->y*mx[1][2] + vec->z*mx[2][2] + vec->w*mx[3][2] ;
  206.     vres->w =
  207.         vec->x*mx[0][3] + vec->y*mx[1][3] + vec->z*mx[2][3] + vec->w*mx[3][3] ;
  208. }
  209.  
  210.  
  211. /*
  212.  * Multiply two 4x4 matrices.
  213.  */
  214. void
  215. lib_matrix_multiply(MATRIX mxres, MATRIX mx1, MATRIX mx2 )
  216. {
  217.     int i, j;
  218.     for ( i = 0; i < 4; i++ ) {
  219.         for ( j = 0; j < 4; j++ ) {
  220.             mxres[i][j] = mx1[i][0]*mx2[0][j] +
  221.                           mx1[i][1]*mx2[1][j] +
  222.                           mx1[i][2]*mx2[2][j] +
  223.                           mx1[i][3]*mx2[3][j] ;
  224.         }
  225.     }
  226. }
  227.  
  228.  
  229. /*
  230.  * Rotate a vector pointing towards the major-axis faces (i.e. the major-axis
  231.  * component of the vector is defined as the largest value) 90 degrees to
  232.  * another cube face.  Mod_face is a face number.
  233.  *
  234.  * If the routine is called six times, with mod_face=0..5, the vector will be
  235.  * rotated to each face of a cube.  Rotations are:
  236.  *     mod_face = 0 mod 3, +Z axis rotate
  237.  *     mod_face = 1 mod 3, +X axis rotate
  238.  *     mod_face = 2 mod 3, -Y axis rotate
  239.  */
  240. void
  241. lib_rotate_cube_face(COORD4 *vec, int major_axis, int mod_face )
  242. {
  243.     double  swap ;
  244.  
  245.  
  246.     mod_face = (mod_face+major_axis) % 3 ;
  247.     if ( mod_face == 0 ) {
  248.         swap   = vec->x ;
  249.         vec->x = -vec->y ;
  250.         vec->y = swap ;
  251.     }
  252.     else if ( mod_face == 1 ) {
  253.         swap   = vec->y ;
  254.         vec->y = -vec->z ;
  255.         vec->z = swap ;
  256.     }
  257.     else {
  258.         swap   = vec->x ;
  259.         vec->x = -vec->z ;
  260.         vec->z = swap ;
  261.     }
  262. }
  263.  
  264.  
  265. /*
  266.  * Portable gaussian random number generator (from "Numerical Recipes", GASDEV)
  267.  * Returns a uniform random deviate between 0.0 and 1.0.  'iseed' must be
  268.  * less than M1 to avoid repetition, and less than (2**31-C1)/A1 [= 300718]
  269.  * to avoid overflow.
  270.  */
  271. #define M1      134456
  272. #define IA1     8121
  273. #define IC1     28411
  274. #define RM1     1.0/M1
  275.  
  276. double
  277. lib_gauss_rand(long iseed)
  278. {
  279.     double  fac ;
  280.     long    ix1, ix2 ;
  281.     double  r ;
  282.     double  v1, v2 ;
  283.  
  284.  
  285.     ix2 = iseed ;
  286.  
  287.     do {
  288.         ix1 = (IC1+ix2*IA1) % M1 ;
  289.         ix2 = (IC1+ix1*IA1) % M1 ;
  290.         v1 = ix1 * 2.0 * RM1 - 1.0 ;
  291.         v2 = ix2 * 2.0 * RM1 - 1.0 ;
  292.         r = v1*v1 + v2*v2 ;
  293.     } while ( r >= 1.0 ) ;
  294.  
  295.     fac = sqrt( (double)( -2.0 * log( (double)r ) / r ) ) ;
  296.     return( v1 * fac ) ;
  297. }
  298.  
  299.  
  300.  
  301.  
  302. /* OUTPUT ROUTINES
  303.  *
  304.  * Files are output as lines of text.  For each entity, the first line
  305.  * defines its type.  The rest of the first line and possibly other lines
  306.  * contain further information about the entity.  Entities include:
  307.  *
  308.  * "v"  - viewing vectors and angles
  309.  * "l"  - positional light location
  310.  * "b"  - background color
  311.  * "f"  - object material properties
  312.  * "c"  - cone or cylinder primitive
  313.  * "s"  - sphere primitive
  314.  * "p"  - polygon primitive
  315.  * "pp" - polygonal patch primitive
  316.  */
  317.  
  318. /*
  319.  * Output viewpoint location.  The parameters are:
  320.  *   From:  the eye location.
  321.  *   At:  a position to be at the center of the image.  A.k.a. "lookat"
  322.  *   Up:  a vector defining which direction is up.
  323.  *
  324.  * Note that no assumptions are made about normalizing the data (e.g. the
  325.  * from-at distance does not have to be 1).  Also, vectors are not
  326.  * required to be perpendicular to each other.
  327.  *
  328.  * For all databases some viewing parameters are always the same:
  329.  *
  330.  *   Viewing angle is defined as from the center of top pixel row to bottom
  331.  *     pixel row and left column to right column.
  332.  *   Yon is "at infinity."
  333.  *   Resolution is always 512 x 512.
  334.  */
  335. void
  336. lib_output_viewpoint( from, at, up, angle, hither, aspect, resx, resy)
  337. COORD4  *from;
  338. COORD4  *at;
  339. COORD4  *up;
  340. double  angle;
  341. double  hither;
  342. double  aspect;
  343. int     resx;
  344. int     resy;
  345. {
  346.     printf( "viewpoint {\n" ) ;
  347.     printf( "   from <%g, %g, %g>\n", from->x, from->y, from->z);
  348.     printf( "   at <%g, %g, %g>\n", at->x, at->y, at->z);
  349.     printf( "   up <%g, %g, %g>\n", up->x, up->y, up->z);
  350.     printf( "   angle %g\n", angle);
  351.     printf( "   hither %g\n", hither);
  352.     printf( "   aspect %g\n", aspect);
  353.     printf( "   resolution %d, %d\n", resx, resy);
  354.     printf( "   }\n");
  355. }
  356.  
  357.  
  358. /*
  359.  * Output light.  A light is defined by position.  All lights have the same
  360.  * intensity.
  361.  *
  362.  */
  363. void
  364. lib_output_light(center_pt)
  365. COORD4  *center_pt ;
  366. {
  367.     printf( "light <%g, %g, %g>\n", center_pt->x, center_pt->y, center_pt->z ) ;
  368. }
  369.  
  370.  
  371. /*
  372.  * Output background color.  A color is simply RGB (monitor dependent, but
  373.  * that's life).  The format is:
  374.  *     "b" red green blue
  375.  */
  376. void
  377. lib_output_background_color( color )
  378. COORD4  *color ;
  379. {
  380.     printf( "background <%g, %g, %g>\n", color->x, color->y, color->z ) ;
  381. }
  382.  
  383. void
  384. lib_output_bounding_slab(COORD4 *dir)
  385. {
  386.     printf("bounding_slab <%g, %g, %g>\n", dir->x, dir->y, dir->z);
  387. }
  388.  
  389.  
  390. /*
  391.  * Output a color and shading parameters for the object in the format:
  392.  *     "f" red green blue Kd Ks Shine T index_of_refraction
  393.  *
  394.  * Kd is the diffuse component, Ks the specular, Shine is the Phong cosine
  395.  * power for highlights, T is transmittance (fraction of light passed per
  396.  * unit).  0 <= Kd <= 1 and 0 <= Ks <= 1, though it is not required that
  397.  * Kd + Ks == 1.
  398.  *
  399.  * The fill color is used to color the objects following it until a new color
  400.  * is assigned or the file ends.
  401.  *
  402.  * Returns the name of the texture that will be used
  403.  */
  404. char *
  405. lib_output_color(color, ka, kd, ks, ang, r, t, i_of_r )
  406. COORD4  *color;
  407. double  ka;
  408. double  kd;
  409. double  ks;
  410. double  ang;
  411. double  r;
  412. double  t;
  413. double  i_of_r;
  414. {
  415.    char *txname = malloc(7*sizeof(char));
  416.    if (txname == NULL) {
  417.       printf("Failed to allocate texture name %d\n", texture_count);
  418.       exit(1);
  419.       }
  420.    sprintf(txname, "txt%03d", texture_count++);
  421.    txname[6] = '\0';
  422.    printf("define %s\n", txname);
  423.    printf("texture {\n");
  424.    printf("   surface {\n");
  425.    printf("      ambient <%g, %g, %g>, %g\n",
  426.           color->x, color->y, color->z, ka);
  427.    printf("      diffuse <%g, %g, %g>, %g\n",
  428.           color->x, color->y, color->z, kd);
  429.    printf("      specular white, %g\n", ks);
  430.    printf("      microfacet Phong %g\n", ang);
  431.    printf("      reflection white, %g\n", r);
  432.    printf("      transmission white, %g, %g\n", t, i_of_r);
  433.    printf("      }\n");
  434.    printf("   }\n");
  435.    return txname;
  436. }
  437.  
  438.  
  439. /*
  440.  * Output cylinder or cone.  A cylinder is defined as having a radius and an
  441.  * axis defined by two points, which also define the top and bottom edge of the
  442.  * cylinder.  A cone is defined similarly, the difference being that the apex
  443.  * and base radii are different.  The apex radius is defined as being smaller
  444.  * than the base radius.  Note that the surface exists without endcaps.
  445.  *
  446.  * If format=OUTPUT_CURVES, output the cylinder/cone in format:
  447.  *     "c"
  448.  *     base.x base.y base.z base_radius
  449.  *     apex.x apex.y apex.z apex_radius
  450.  *
  451.  */
  452. void
  453. lib_output_cylcone( base_pt, apex_pt, texture_name)
  454.    COORD4 *base_pt, *apex_pt;
  455.    char *texture_name;
  456. {
  457.    printf("object {\n");
  458.    if (base_pt->w == apex_pt->w)
  459.       printf("   cylinder <%g, %g, %g>, <%g, %g, %g>, %g\n",
  460.              base_pt->x, base_pt->y, base_pt->z,
  461.              apex_pt->x, apex_pt->y, apex_pt->z, apex_pt->w);
  462.    else
  463.       printf("   cone <%g, %g, %g>, %g, <%g, %g, %g>, %g\n",
  464.              base_pt->x, base_pt->y, base_pt->z, base_pt->w,
  465.              apex_pt->x, apex_pt->y, apex_pt->z, apex_pt->w);
  466.    printf("   %s\n", texture_name);
  467.    printf("   }\n");
  468. }
  469.  
  470.  
  471. /*
  472.  * Output sphere.  A sphere is defined by a radius and center position.
  473.  */
  474. void
  475. lib_output_sphere( center, texture_name )
  476.    COORD4 *center;
  477.    char *texture_name;
  478. {
  479.    printf("object {\n");
  480.    printf("   sphere <%g, %g, %g>, %g\n",
  481.           center->x, center->y, center->z, center->w);
  482.    printf("   %s\n", texture_name);
  483.    printf("   }\n");
  484. }
  485.  
  486.  
  487. /*
  488.  * Output polygon.  A polygon is defined by a set of vertices.  With these
  489.  * databases, a polygon is defined to have all points coplanar.  A polygon has
  490.  * only one side, with the order of the vertices being counterclockwise as you
  491.  * face the polygon (right-handed coordinate system).
  492.  *
  493.  * The output format is always:
  494.  *     "p" total_vertices
  495.  *     vert1.x vert1.y vert1.z
  496.  *     [etc. for total_vertices polygons]
  497.  *
  498.  */
  499. void
  500. lib_output_polygon(tot_vert, vert, texture_name)
  501.    int tot_vert;
  502.    COORD4 *vert;
  503.    char *texture_name;
  504. {
  505.    int num_vert;
  506.    printf("object {\n");
  507.    printf("   polygon %d,\n", tot_vert);
  508.    for (num_vert = 0; num_vert < tot_vert; num_vert++) {
  509.       printf("      <%g, %g, %g>",
  510.              vert[num_vert].x, vert[num_vert].y, vert[num_vert].z);
  511.       if (num_vert < tot_vert-1)
  512.          printf(",");
  513.       printf("\n");
  514.       }
  515.    printf("   %s\n", texture_name);
  516.    printf("   }\n");
  517. }
  518.  
  519.  
  520. /*
  521.  * Output polygonal patch.  A patch is defined by a set of vertices and their
  522.  * normals.  With these databases, a patch is defined to have all points
  523.  * coplanar.  A patch has only one side, with the order of the vertices being
  524.  * counterclockwise as you face the patch (right-handed coordinate system).
  525.  *
  526.  * The output format is always:
  527.  *     "pp" total_vertices
  528.  *     vert1.x vert1.y vert1.z norm1.x norm1.y norm1.z
  529.  *     [etc. for total_vertices polygonal patches]
  530.  *
  531.  */
  532. void
  533. lib_output_polypatch(tot_vert, vert, norm, texture_name)
  534.    int tot_vert;
  535.    COORD4 *vert;
  536.    COORD4 *norm;
  537.    char *texture_name;
  538. {
  539.    if (tot_vert != 3) {
  540.       printf("Patches must have 3 vertices, input was %d\n", tot_vert);
  541.       exit(1);
  542.       }
  543.    printf("object {\n");
  544.    printf("   patch <%g, %g, %g>, <%g, %g, %g>, <%g, %g, %g>,\n",
  545.           vert[0].x, vert[0].y, vert[0].z,
  546.           norm[0].x, norm[0].y, norm[0].z,
  547.           vert[1].x, vert[1].y, vert[1].z);
  548.    printf("         <%g, %g, %g>, <%g, %g, %g>, <%g, %g, %g>\n",
  549.           norm[1].x, norm[1].y, norm[1].z,
  550.           vert[2].x, vert[2].y, vert[2].z,
  551.           norm[2].x, norm[2].y, norm[2].z);
  552.    printf("   %s\n", texture_name);
  553.    printf("   }\n");
  554. }
  555.